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Abstract 

The majority of practical flows, particularly those flows in applications 
of importance to transport, distribution and climate, are turbulent and as 
a result experience complex three-dimensional motion with increased drag 
compared with the smoother, laminar condition. In this study, we describe 
the development of a simple model that predicts important structural and 
scaling features of wall turbulence. We show that a simple linear superposi- 
tion of modes derived from a forcing-response analysis of the Navier-Stokes 
equations can be used to reconcile certain key statistical and structural de- 
scriptions of wall turbulence. The computationally cheap approach explains 
and predicts vortical structures and velocity statistics of turbulent flows that 
have previously been identified only in experiments or by direct numerical 
simulation. In particular, we propose an economical explanation for the 
meandering appearance of very large scale motions observed in turbulent 
pipe flow, and likewise demonstrate that hairpin vortices are predicted by the 
model. This new capability has clear implications for modeling, simulation 
and control of a ubiquitous class of wall flows. 

In flows of the most widely-ranging practical interest, namely flows over sur- 
faces or wall turbulence, the turbulence problem is greatly complicated by inho- 
mogeneity in the wall-normal direction caused by the no-slip and no-penetration 
boundary conditions at the wall. As such, there has been an enduring practical and 
intellectual interest in understanding how to describe the state of a range of turbu- 
lent flows with a view to controlling and optimizing their behavior, which could 
help reduce drag losses in applications like transport vehicles or pipelines. 

The nonlinear equations of motion for fluids were written down by Navier and 
Stokes in the early 19th century. However turbulence, a phenomenon involving 
interaction at a very broad range of spatial and time scales, has proved difficult 
to model without resorting to expensive direct numerical simulation. Turbulent 
fluctuations of velocity and pressure have characteristics of chaotic motion, and 
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a common feature of wall turbulence is the presence of persistent, coherent vor- 
tical structures which inhabit distinct regions of the flow ll33l l32l and have been 
generally believed to be dominated by nonlinear dynamics. These observations of 
distinct classes of coherent structure in parallel with the measurement of turbu- 
lent velocity statistics over the last sixty or more years have led to a dichotomy in 
the understanding of the formation, development and scaling of even the simplest, 
canonical turbulent wall flows, such as flow through straight channels and pipes, or 
flow over a flat plate in the absence of a pressure gradient. 

In this work, we demonstrate that the formation of such coherent structures is a 
natural consequence of the velocity field that arises as a near-singular response of 
the linearised equations of motion to the wave-like forcing arising due to nonlin- 
ear interactions with other wave-like motions. The framework provides predictive 
information both about how the wall-normal distribution of turbulent energy of the 
velocity field shows self- similarity across Reynolds numbers, and also describes 
the wall-normal coherence of structures within the flow. The predictions of the 
scaling and distribution of second-order velocity statistics, and of structural char- 
acteristics are provided using an essentially linear model. The analysis uses only 
the Navier-Stokes equations and an assumed mean velocity profile. 

The existence of coherent vortical structures in wall turbulence has been known 
for many decades [33], but a quantitative or ab initio mathematical means of de- 
scribing the full flow field in terms of these building blocks has been lacking. There 
have been many observations of a common structural pattern described as a hair- 
pin vortex (a vortical loop with legs originating close to the wall, a body inclined 
in the downstream direction and a sense of rotation consistent with the vorticity 
associated with the mean shear) in both experiment and simulation |fT3ll3ll2l[36ll, 
while the attached eddy hypothesis formulated by Townsend [34] and Perry and 
co-authors ||2"9l in the mid- 20th century has sought to use phenomenological ar- 
guments to determine the velocity field associated with hierarchies of these struc- 
tures. Nonlinear mechanisms for the formation and growth of packets of hairpins 
in otherwise quiescent flow have also been proposed 1381 . while there remains con- 
troversy as to whether they are simply a remnant of the transition to turbulence, or 
even whether whole hairpins, rather than a statistical vortical imprint, even exist. 

With alternative models lacking, statistical descriptions have remained the sim- 
plest method to obtain comprehensive descriptions of the turbulent field. Under the 
assumptions of stationarity and ergodicity, simple spectral representations provide 
information on important scales in the flow. While information at streamwise and 
spanwise wavenumber and temporal frequency (k,n, CO, respectively) is required to 
fully describe the flow, usually one-dimensional (integrated) spectra are reported 
due to the intensive nature of obtaining and storing information in three dimen- 
sions. The streamwise spectra reveal the presence of the widest range of scales, 
from the Kolmogorov scales |[T8l[T9ll responsible for small-scale dissipation of en- 
ergy up to so-called very large scale motions (VLSMs) of the order of ten times 
the length scale imposed by the flow geometry (typically the pipe radius, or the 
boundary layer thickness). Some researchers have inferred that the VLSMs reach 
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lengths of the order of thirty radii when spanwise meander of the coherent region 
is taken into account EoTl . a topic of continuing interest in the literature. A non- 
dimensional representation of ratio of the scales that are important far from the wall 
and in the near-wall region is given by the Reynolds number Re, a parameter that 
rapidly increases with increasing velocity^] A recent focus of research has been 
on the origin and development of large and small scales, and how their competing 
influences lead to trends with increasing Reynolds number. While the origin of the 
VLSMs has remained elusive, there has been recent progress in predicting velocity 
statistics by considering the nature of their spectral interaction with the energetic 
near- wall region [20]. 

A truly complete interpretation of a turbulent flow field requires assimilation of 
both the velocity statistics and vortical structures; observational progress has been 
made by considering the generation of one by the other, posing a classic "chicken 
and egg" conundrum. In the absence of low-order or simple predictive models for 
turbulent flows, the only current recourse for the fluid mechanician who wishes to 
quantify the state of turbulence lies in expensive simulations or experiments. 

One recent and promising approach to understanding the development of the 
turbulent field has been the analysis of the receptivity of wall flows to particular 
forcing or disturbances. In the case of transition to turbulence from laminar flow, 
progress has been made PTTl by considering the very large amplification that is 
possible arising from the evolution of finite perturbations to the laminar flows. 
This approach is best suited to explaining sub-critical transition to turbulence and 
it is unclear how suitable such an analysis is in fully turbulent flow. The approach 
generalises straightforwardly to complex geometries [1] and the current approach 
should also generalise in such a way. 

There has been some investigation of the non-normal growth mechanisms in 
turbulent flow, for instance the early attempt by @ to predict the spacing of near- 
wall streaks in turbulent flow and the more recent studies of [10], [9] and |[35l . 
These approaches exploit the properties of the linear operator to "select" by prefer- 
ential amplification the "optimal" disturbance that leads to the largest amplification 
in some norm across all disturbances. A notable and somewhat limiting issue in 
the treatment of turbulent flow is the modelling of the interaction of the amplified 
disturbances with the "background" turbulence. One solution proposed OUl was 
to use the eddy viscosity formulation of Q ; this relies on an a priori knowledge 
of the spatially-averaged, wall-normal variation of the mean Reynolds stress inte- 
grated across contributions from various Reynolds numbers. Other attempts have 
been made to use linear analysis to explain the dominant features of turbulent flow 
in terms of optimal transient modes in the initial value problem B[T0j[9l, response 
to stochastic forcing lfl2l l4il and system norm analysis lfl6l l25l . In recent work, 

1 The Reynolds number described here is defined by the ratio between a flow length scale (here 
the pipe radius R) and a scale representative of small, near-wall motions, the viscous length scale 
defined by V ju x , where v is the kinematic viscosity of the fluid and u x is a velocity scale associated 
with the skin friction acting on the wall, T, and the fluid density. Thus our Reynolds number is given 
by Re = R + = Ru t /v (alternative definitions are possible). 
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|[35l have investigated the maximal response to harmonic forcing in pipe flow. Per- 
haps most importantly, it has been shown that both linear non-normality lfl4l and 
the terms that are linear in the turbulent fluctuation ifTTl are required to sustain 
turbulence in infinite or periodic wall-bounded flows. 

Alternatively it is possible to examine the response of the linearised system to 
forcing. In this picture, the observed behaviour is explained by even small forc- 
ing on the system leading to energetic flowfield response (41 IT6J . The device of 
identifying the nonlinear interaction between Fourier modes in the Navier-Stokes 
equations (NSE) as a forcing actually acting on the a linear system permits the 
extension of such methods to fully developed turbulent flows [24], allowing the 
successful prediction of observed features of turbulent flow. 

We consider turbulent flow through a long straight pipe with a cylindrical cross- 
section. Laminar flow in this geometry is stable to infinitesimal disturbances, and 
the transition to turbulence is still not completely understood [15], but the pipe 
offers the analytical benefits of statistical homogeneity in the streamwise direction 
and a simple constraint on the azimuthal wavenumber. Turbulent flow through 
pipes is important for applications like the transport of fluids such as oil and natural 
gas, and also in numerous natural and biomedical applications, and is also highly 
relevant to the study of other canonical flows. 

We recently formulated a traveling wave framework by which to analyze the 
dominant velocity mode shapes at particular wavenumber-frequency combinations 
ll24l . The fully turbulent velocity field, v, can be represented in a divergence- 
free basis as a superposition of Fourier modes at various spatial wavenumbers and 
temporal frequencies with wall-normal variation (waves), 



where y is the wall-normal distance, k and n are the wavenumbers in the stream- 
wise (x) and azimuthal (8) directions, all normalized with the pipe radius, and 
CO is the temporal frequency with respect to time t. These waves, which are he- 
lical in pipe flow, represent the transmission of fluctuating energy relative to the 
mean flow. Under this decomposition, the Navier-Stokes equations can be writ- 
ten in an input-output formulation, where the characteristics of the linear operator 
JZk,n,ca(y) describe the velocity field's response to an unmodelled harmonic forcing 
f, itself arising from the nonlinear interactions between the velocity field at other 
wavenumber-frequency combinations. The equations expressed in this form are 



and the operator (ico — Jj?k,n,a(y))~ l relating the forcing to the velocity field re- 
sponse is called the resolvent. 

The analysis up to this point is very similar to the development of the linearized, 
fourth-order Orr-Sommerfeld-Squire operator of linear stability theory, with the 
exceptions that the (k,n,co) = (0,0,0) mode is identified as the turbulent mean 




\(y,k,n,co) = (ico - ^k^caiy)) l f(y,k,n,co) 
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velocity rather than the laminar solution, and the nonlinear forcing terms are ex- 
plicitly retained in the present analysis. As such, receptivity concepts relevant to 
the study of neutrally-stable disturbances in inviscid, linearized laminar flow can 
be extended to the turbulent case, with the understanding that in the latter case the 
waves are lightly damped and would asymptotically decay in the absence of forc- 
ing f. The form of the resolvent dictates that a nearly-singular, essentially inviscid 
response occurs at the critical layer, where the local mean velocity U is equal to 
the streamwise convective velocity of the wave. 

Given the very selective receptivity of the flow, we may hypothesise that the 
velocity field is dominated by the largest possible response to forcing, at any given 
wavenumber and frequency set. Mathematically, we may find and order these or- 
thogonal response modes by use of the Schmidt decomposition ||37ll of the resol- 
vent, which provides, for a given forcing magnitude, the largest possible velocity 
field response, the second largest and so onj^] In practice, we find that the velocity 
field given by the first of these response modes is associated with a response typ- 
ically one to many orders of magnitude larger than that of the second mode. The 
velocity field observed in a real flow will therefore be well described by the first 
response mode if the forcing contains a non-negligible component of the correct 
shape. This assumption of selective receptivity is congruent with an assumption 
that the resolvent may be approximated by a low rank operator (which has been 
postulated in various forms by many previous researchers). 

Under the further assumption of forcing that is small relative to the magnitude 
of the response modes, a standard asymptotic analysis tfTTl can be used to describe 
the scaling of the two regions where viscous effects are required: firstly at the 
critical layer (to regularize the singularity of the inviscid problem) and at the wall 
(to meet the wall boundary condition). We call a response mode in which viscous 
modification at the wall dominates, a wall or attached mode, and we call a response 
mode where viscous modifications at the critical layer dominate a critical mode. 

In contrast to calculating statistics at a single wall-normal location (y) at a given 
Fourier wavenumber, the response modes provide a basis that predicts the wall- 
normal amplitude and phase variation of the velocity field. To date, we have exam- 
ined the structure of the response modes up to a Reynolds number of Re = 2 x 10 , 
of the same order of magnitude as conditions in, for example, a transcontinental 
natural gas pipeline [22]. Here we present key results at Re = 2 x 10 3 , a condition 
that is achievable in both experiment and state-of-the-art simulation. 

Figure[T]shows the streamwise velocity component arising from a combination 
of the left and right going helical velocity response modes with a wavenumber- 
frequency combination representative of the dominant motion near the wall, namely 
streamwise and spanwise wavelengths of a thousand and one hundred viscous 
units, respectively, and a convection velocity of ten times the friction velocity, 
A z + , t/+) = (1000, 100, 10). The distinctive pattern of rolling motions aligned 

2 The calculations are performed using a modified version of the approach of |25 1 and the experi- 
mental mean velocity data of [23|. 
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Figure 1 : Shape of the first singular mode representative of the dominant near wall 
motions, (A v + , {/+) = (1000, 100,10). Color denotes isosurfaces of streamwise 
velocity (streaks), where red and blue correspond to high and low velocity respec- 
tively relative to the mean flow (heading into the page), and the arrows show the 
sense of the in-plane velocity field. 

in the streamwise direction and strong, alternating inclined streaks of fast and slow 
streamwise velocity u shown in Figure[T]is entirely consistent with visual and quan- 
titative observations of the near-wall region in canonical flows [32]. In McKeon & 
Sharma ll24l . we showed that the wall- normal location of the peak intensity of u 
associated with this mode is independent of Reynolds number. This scaling result 
is borne out by experimental measurements over a range of Reynolds numbers. 

Our framework can also be used predictively to describe the Reynolds number 
dependence of the location of peak VLSM energy, yi, and to understand the origin 
of VLSMs. A mode that is both attached to the wall and critical has a special 
significance. In pipe flow, for parameters representative of a VLSM, (k,n) = (l, 10) 
approximately, this condition occurs when the convective velocity is 2/3 of the 
centerline velocity, independent of Reynolds number, leading to the prediction 

y+ =0.%Re 2 l\ 

The agreement between this theoretically derived expression and the experimental 
results of E71 reported by McKeon [21] is remarkably good II241 . This prediction 
differs from the Re 1 ' 2 law appropriate for a boundary layer. 

Having described important statistical aspects of pipe flow, we now turn our at- 
tention to structural considerations. A simple superposition of first response modes 
at different wavenumber-frequency combinations also sheds light on the contro- 
versy surrounding the true length of the VLSM motions. Even the addition of the 
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Figure 2: Regeneration of apparent "meandering" of the VLSMs. The top panel 
contains isosurfaces of the streamwise velocity at y/R = 0.15 for the summa- 
tion of left- and right-going VLSM modes and the middle panel shows how even 
longer coherence, of the order of the panel length, can be obtained by superposing 
the velocity response mode pair fields from two additional (shorter) modes with 
(k,n,a>) = (4.7, ±12,0.2) and (6.2, ±15, 0.6) with amplitude 75% of the VLSM 
modes. The bottom panel shows the experimental results of Monty et al. [26). 
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velocity fields corresponding to only two additional pairs of response modes with 
similar amplitudes to the VLSM mode described above quickly leads to the obser- 
vation of apparently meandering structures with length far greater than six radii, 
as shown in Figure [2] The visual similarity is striking to experimental data shown 
in the lowest panel of the figure (from GoTl ). The meandering phenomenon is re- 
vealed to be an artefact of the many response modes that are present in a real flow 
combining with the energetic content of the VLSM itself, effectively decorrelating 
the VLSM mode. 

The three-dimensional velocity field associated with the attached near-wall 
mode shown in Figure [T] gives an intuitive hint as to the locations of coherent 
vorticity associated with this type of mode. We identify structure through mea- 
sures which distinguish between the shear and rotational components of vorticity, 
namely the symmetric and anti-symmetric components of the velocity gradient ten- 
sor, Vv. While any of the commonly used measures (see Chakraborty et al. (H) 
give very similar results, we highlight isosurfaces of swirling strength^for a wall 
mode with (k,n,(o) = (4.5, ±10, 1.67) in Figure [3^. The surfaces are colored by 
the magnitude of the azimuthal vorticity, with blue corresponding to prograde vor- 
tices, with the same sense of rotation as the traditional hairpin, and red to retro- 
grade vortices, with the opposite sense of rotation. The latter have been reported in 
the literature to occur relatively infrequently I28l l6l. Our model predicts an equal 
distribution of prograde and retrograde vortices associated with each wall mode. 
However, in the presence of a mean velocity profile with decreasing shear in the 
wall-normal direction, the retrograde vortices are suppressed while the prograde 
ones are reinforced, leading to a distribution of vortices that is consistent with ex- 
perimental observations, Figure [3}). In this figure, the mean shear is sufficiently 
strong to completely suppress the retrograde vortices below the swirl threshold se- 
lected for plotting. We emphasize that the swirl field regenerated from the full 
range of modes projected out, for example, from a direct numerical simulation 
would reproduce the full swirl field: we explore here a decomposition of the swirl 
field around the mean shear, which precludes consideration of the contribution of 
hairpin heads to the mean shear itself, as explored by other authors, e.g. 0- There 
are two points to note here, firstly, this is a simple manifestation of the mean shear 
being the only net source of azimuthal vorticity, and secondly, swirling strength 
is not a distributive operation. In other words, swirl(a + b) / swirl(a) + swirl(ft). 
Thus the phenomenon is a direct consequence of the diagnostics commonly used 
to identify structure. To illustrate, consider the example of swirling strength in 
a two-dimensional, streamwise/wall-normal velocity field (u,v). Examination of 
the expression for X c i given in the equation below shows that the mean shear will 
lead to a reduction in apparent swirling strength in regions where dv/dx > 0, as 

3 A C ,-, or the imaginary part of the complex conjugate eigenvalue pair associated with the velocity 
gradient tensor 
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Figure 3: Isosurfaces of constant swirling strength for the (k,n,<a) = 
(4.5, ±10, 1.67) velocity response mode (three wavelengths are shown in the 
streamwise and two in the spanwise directions) at Re = 1800, color-coded with 
the sense of the azimuthal rotation. Blue and red denote pro- and retro-grade swirl 
(or rotation in and counter to the sense of^the classical hairpin vortex), respectively. 
Top: under a Galilean transformation (i.e. constant convection velocity subtracted 
throughout the field of view) there are even numbers of prograde and retrograde 
vortices. Bottom: with the mean velocity profile added, the retrograde vortices 
disappear and the prograde ones are strengthened. 




Figure 4: Experimental streamwise and wall-normal velocity fields obtained us- 
ing PIV in a zero pressure gradient turbulent boundary layer at a similar Reynolds 
number to the results of Figure [3] Top: Galilean transformation (i.e. constant 
convection velocity subtracted throughout the field of view). Bottom: Reynolds 
decomposition (i.e. local mean velocity subtracted at each wall-normal location). 
Blue and red denote pro- and retro-grade swirling strength, respectively, that ex- 
ceeds the chosen threshold magnitude and the swirl is overlaid on a quiver plot of 
the local in-plane velocities. 
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associated with retrograde vortices. 



^ lj // du dv\ 2 ^fdudv dudv 
" 2 y \dx dy J \dx dy dy dx 

Experimental velocity fields obtained in a zero pressure gradient turbulent bound- 
ary layer using particle image velocimetry (PIV) confirm this observation. Figure 
[4^ shows the more usual Galilean decomposition of the data (where a constant 
convection velocity has been subtracted from all locations in the field). There is a 
prevalence of prograde vortices, denoted in blue, whereas the Reynolds-decomposed 
field of Figure [4j>, in which the local mean velocity has been subtracted from all 
wall-normal locations, shows an essentially equal number of prograde and retro- 
grade cores. The positive local du/dy associated with the most energetic mode 
in each picture explains the apparent alignment of hairpins into hairpin packets 
identified with ramp-like shear layers, while the relative motion between fast- and 
slow-moving modes explains the apparent wall-normal growth and decay of indi- 
vidual vortices, and renders the concept of relative phase between modes a moot 
point. The linearity of our model reveals the crucial importance of a swirl diagnos- 
tic that is not affected by mean shear, a failing that has impeded our understanding 
of the appropriate distribution of coherent vortical structures that corresponds to 
the appropriate mathematical understanding. 

It is clear that a very complex velocity field can be obtained simply by super- 
posing modes with different (k,n, co) (and therefore different convective velocities) 
and amplitudes, and with associated local velocity gradients. Critically, because 
swirling strength is not a distributive operation, the swirl field cannot be simply 
determined from a superposition of the swirling strength associated with the indi- 
vidual response modes. 

Understanding turbulence has proven difficult because it requires modeling at 
all scales. Our approach predicts and explains aspects of the vortical structure 
and velocity statistics associated with turbulent flows that previously have been 
explained only in phenomenological terms and have been identified only in exper- 
iments or computationally costly direct numerical simulation, which is limited to 
Reynold numbers much lower that those relevant to most applications. The linear- 
ity of the processes driving wall turbulence permits superposition of the response 
modes, a surprising result given the understanding of the importance of nonlin- 
ear dynamics in wall turbulence, although that linear processes are also important 
is well known. This linearity is revealed by formulating the NSE as a forcing- 
response problem. We have studied pipe flow in the current work, but the same 
approach can be applied to both internal and external flows with simple modifi- 
cations. The results shown in this paper were generated in seconds using a stan- 
dard laptop computer, and the approach extends to higher Reynolds numbers. The 
only limitation is the numerical precision required to deal with the near-singular 
response to the most amplified forcing. 

A fundamental understanding of the underlying structure of wall turbulence 
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has eluded researchers for decades. Excitingly, our results close in on a recon- 
ciliation of the statistical and structural interpretations of such flows by working 
from the NSE and an assumed mean profile. The understanding of the different 
types of mode and where they occur in the flow has significant implications for 
the prediction of Reynolds number trends and the modeling of turbulent activity 
at reduced computational cost. While the forcing in our model is currently un- 
structured, such that we do not determine the appropriate amplitudes at individual 
modes, the possibility to "close the loop" and formulate a reduced order model of 
turbulent pipe flow is compelling: under a self-consistent combination of modes 
the assumed mean velocity profile will be generated, such that a self-sustaining 
system can be designed. The potential to design rigorous control techniques for 
such a system with the objective of enforcing favorable turbulence characteristics 
is a natural, and plausible, next step. 
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